{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/google/earthengine-api/blob/master/demos/flooded-roads/ExportToBigQuery.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "wE93poeBnMnW"
      },
      "source": [
        "# Demo Scenario\n",
        "\n",
        " Extreme weather events have a devastating impact around the world. Flooding, heat waves, and drought have substantial human and financial costs, causing mortality and devastation of homes and property. The following example shows how to use satellite data mosaics from Earth Engine and open road datasets from BigQuery, processing the data in both environments to determine which road segments are affected by a flooding event in the UK.\n",
        "\n",
        "\n",
        " ## Demo Flow\n",
        "\n",
        " Earth Engine => Big Query => Visualization\n",
        "\n",
        " ## Prerequisites\n",
        "\n",
        "\n",
        "1.   Install the dependencies {--geemap , earthengine-api}\n",
        "2.   Create a new GCP project and check that billing is enabled. Alternatively an existing project with valid billing ID could be used.\n",
        "3.   Enable the BigQuery and Earth Engine APIs.\n",
        "4.   Make sure to have the permissions to create dataset.\n",
        "5.   Make sure to have the below IAM roles assigned on the BigQuery dataset\n",
        "\n",
        ">1. bigquery.dataEditor + bigquery.jobUser\n",
        "\n",
        ">2. bigquery.dataOwner + bigquery.jobUser\n",
        "\n",
        ">3. bigquery.admin\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "352zE5rZDO0v"
      },
      "outputs": [],
      "source": [
        "# @title Prerequisites\n",
        "# Install 'geemap' library to display the map.\n",
        "!pip install geemap\n",
        "!pip install earthengine-api --upgrade"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "llYex6u5jNGx"
      },
      "outputs": [],
      "source": [
        "# @title Parameter Setup\n",
        "# Replace the project id with your project.\n",
        "project_id = \"example-project\"  # @param {type:\"string\"}\n",
        "dataset_id = \"ee_export\"\n",
        "table_id = \"ee_test\"\n",
        "table =  dataset_id + \".\" + table_id\n",
        "region = 'us'\n",
        "table_path = project_id + \".\" + dataset_id + \".\" + table_id\n",
        "print(\"Region:      \",region)\n",
        "print(\"Table Path:  \",table_path)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "MQJDgFB0jEH0"
      },
      "outputs": [],
      "source": [
        "# @title Authenticate and initialize the Session\n",
        "import google\n",
        "import ee\n",
        "from google.cloud import bigquery\n",
        "from google.colab import auth as google_auth\n",
        "client= bigquery.Client()\n",
        "google_auth.authenticate_user()\n",
        "credentials, auth_project_id = google.auth.default()\n",
        "ee.Initialize(credentials, project=project_id)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "9guXwDgP-grd"
      },
      "outputs": [],
      "source": [
        "# @title Create BQ Dataset\n",
        "# Create a BQ dataset.\n",
        "!bq --location={region} mk --dataset {project_id}:{dataset_id}"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "nNcSsY9y5TAK"
      },
      "source": [
        "Define area of interest and centre point"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "o2QgpHD5-PRG"
      },
      "outputs": [],
      "source": [
        "# @title Define point and Area of Interest(aoi)\n",
        "# Define AOI (Area of Interest) polygon\n",
        "aoi = ee.Geometry.Polygon([[-2.92, 54.10],\n",
        "          [-2.92, 53.99],\n",
        "          [-2.67, 53.99],\n",
        "          [-2.67, 54.10]])"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "8wIR_Fwu04En"
      },
      "outputs": [],
      "source": [
        "# @title Display AOI and point\n",
        "import geemap\n",
        "Map = geemap.Map()\n",
        "Map.centerObject(aoi, 12);\n",
        "Map.setOptions(mapTypeId='HYBRID', styles={}, types=[])\n",
        "Map.addLayer(aoi, {\"color\":\"blue\"});\n",
        "Map"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ycom7Zl36R4I"
      },
      "source": [
        "The Earth Engine Data Catalog contains the [Copernicus Sentinel Synthetic Aperture Radar collection](https://colab.research.google.com/drive/11Cr9nqizvw4D-SwSKgV2mZ6njFiy1lac#scrollTo=ycom7Zl36R4I&line=1&uniqifier=1). This public dataset is composed of radar images that measure how surfaces scatter light waves back to a satellite's sensor. Standing bodies of water act like mirrors for radio signals, reflecting the satellite's radar light away rather than scattering it back to the imaging sensor. Most natural surfaces don't have this property, which means that one can differentiate standing bodies of water from their surroundings by looking for \"dark\" patches in the images (that is, areas with low backscatter values). Let’s prepare the input data by selecting an area of interest and filtering images with vertical-vertical (\"VV\") polarization, sending vertically polarized light, and measuring the vertically polarized light that's returned."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "43gLFRiwKjWg"
      },
      "outputs": [],
      "source": [
        "# @title Data Collections\n",
        "# Load Sentinel-1 C-band SAR Ground Range collection (log scaling, VV co-polar).\n",
        "collection = ee.ImageCollection('COPERNICUS/S1_GRD') \\\n",
        "    .filterBounds(aoi) \\\n",
        "    .filter(ee.Filter.listContains('transmitterReceiverPolarisation', 'VV')) \\\n",
        "    .select('VV')\n",
        "# Smooth the data to remove noise.\n",
        "smoothing_radius = 100  # meters\n",
        "# Filter by date.\n",
        "before = collection.filterDate('2017-11-01', '2017-11-17') \\\n",
        "    .mosaic() \\\n",
        "    .focal_median(smoothing_radius, 'circle', 'meters')  # before floods\n",
        "after = collection.filterDate('2017-11-18', '2017-11-23') \\\n",
        "    .mosaic() \\\n",
        "    .focal_median(smoothing_radius, 'circle', 'meters')  # after floods"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "5hXZvoMYKmbK"
      },
      "outputs": [],
      "source": [
        "# @title Identify Flooded Areas\n",
        "# Threshold smoothed radar intensities to identify areas with standing water.\n",
        "diff_upper_threshold = -3  # dB\n",
        "diff_smoothed = after.subtract(before);\n",
        "diff_thresholded = diff_smoothed.lt(diff_upper_threshold)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "4HCLn7O1nUiQ"
      },
      "source": [
        "use the [Global Surface Water dataset](https://developers.google.com/earth-engine/tutorials/tutorial_global_surface_water_01) to remove persistent surface water (like lakes, rivers, etc.) from the result:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "EzTNlnZtTIEY"
      },
      "outputs": [],
      "source": [
        "# @title Remove Water Areas(other than floods)\n",
        "# Remove global surface water (oceans, lakes, etc.).\n",
        "jrc_data0 = ee.Image(\"JRC/GSW1_0/Metadata\") \\\n",
        "    .select('total_obs') \\\n",
        "    .lte(0)\n",
        "water_occ = ee.Image(\"JRC/GSW1_0/GlobalSurfaceWater\") \\\n",
        "    .select('occurrence') \\\n",
        "    .unmask(0) \\\n",
        "    .max(jrc_data0) \\\n",
        "    .lt(10)\n",
        "diff_thresholded = diff_thresholded.updateMask(water_occ)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "vkowcJLHIZTm"
      },
      "outputs": [],
      "source": [
        "# @title Visualize the Map with Flooded Area\n",
        "# Display flooded areas on the map.\n",
        "import geemap\n",
        "vis_params = {\n",
        "    \"palette\": [\"blue\"],\n",
        "}\n",
        "Map = geemap.Map()\n",
        "Map.setOptions(mapTypeId='HYBRID', styles={}, types=[])\n",
        "Map.centerObject(aoi, 12);\n",
        "Map.addLayer(\n",
        "    diff_thresholded.updateMask(diff_thresholded),\n",
        "    vis_params,\n",
        "    'flooded areas - blue',\n",
        "    True)\n",
        "Map"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "qCJ1tiPG7SBd"
      },
      "source": [
        "We want the flooded areas in BigQuery, so let’s convert flooded pixel data to vector format."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "bwooSOpwSWne"
      },
      "outputs": [],
      "source": [
        "# @title Extract Vectors from the Flooded areas\n",
        "# Extract vectors from the diff threshold to load to BigQuery.\n",
        "vectors = diff_thresholded.reduceToVectors(\n",
        "    geometry = aoi,\n",
        "    scale = 10,\n",
        "    geometryType = 'polygon',\n",
        "    eightConnected = False)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "lrajMu4P8MYr"
      },
      "source": [
        "This is where our new BigQuery connector simplifies export to just making one call: Export.table.toBigQuery."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "cB4s6DxZNqTn"
      },
      "outputs": [],
      "source": [
        "# @title Export to BigQuery\n",
        "task_config = {\n",
        "  'collection': vectors,\n",
        "  'description':'ee2bq_export_polygons',\n",
        "  'table': table_path\n",
        "}\n",
        "task = ee.batch.Export.table.toBigQuery(**task_config)\n",
        "task.start()\n",
        "# The task should run for about a minute. Check task.status() to see the result.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "J2quoziwTDnN"
      },
      "outputs": [],
      "source": [
        "# @title Check Export Job Status\n",
        "# Check the results and make sure the status is COMPLETED before checking the\n",
        "# results in BigQuery.\n",
        "task.status()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ivAhkbOT8wuQ"
      },
      "source": [
        "Now we have exported data available in BigQuery. Execute the query to check the data"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "38K0qKgKxgBI"
      },
      "outputs": [],
      "source": [
        "# @title Check BigQuery Results\n",
        "%%bigquery --project $project_id\n",
        "SELECT * from ee_export.ee_test\n",
        "# If you get a \"table not found\" error then check the step above to see if your\n",
        "# job has completed."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "GhGpv47M9N6Y"
      },
      "source": [
        "In our example we are going to use the public “planet_ways” dataset from OpenStreetMap, so we can find roads that are underwater.\n",
        "\n",
        "\n",
        "> Get polygons corresponding to flat areas from the exported table.\n",
        "\n",
        "> Filter out administrative areas\n",
        "\n",
        "> Join result with road polygons from OpenStreetMap data that intersect with flat areas.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "030QbjYC9M3L"
      },
      "outputs": [],
      "source": [
        "# @title BQ Result set to display flooded highways\n",
        "%%bigquery regions_by_country --project $project_id\n",
        "SELECT\n",
        " id, area,version,changeset,osm_timestamp,ST_ASGEOJSON(flood_poly) as flood_poly,\n",
        " ST_ASGEOJSON(road_geometry) as road_geometry\n",
        "FROM (\n",
        " -- query 1 - find all the flooding areas\n",
        " SELECT\n",
        "   geo AS flood_poly,\n",
        "   ST_AREA(geo) AS area\n",
        " FROM\n",
        "   ee_export.ee_test\n",
        " WHERE\n",
        "   ST_AREA(geo) < 500000 ) t1 -- eliminate admin areas in the dataset\n",
        "JOIN (\n",
        " SELECT\n",
        "   id,\n",
        "   version,\n",
        "   changeset,\n",
        "   osm_timestamp,\n",
        "   geometry as road_geometry\n",
        " FROM\n",
        "   `bigquery-public-data.geo_openstreetmap.planet_ways` planet_ways,\n",
        "   planet_ways.all_tags AS all_tags\n",
        " WHERE\n",
        "   all_tags.key = 'highway' )\n",
        "ON\n",
        " ST_INTERSECTS(flood_poly, road_geometry)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "SjzGS6pF_yqL"
      },
      "outputs": [],
      "source": [
        "# @title Extract the Features from flood_poly\n",
        "floods = '{\"type\": \"FeatureCollection\", \"features\":['\n",
        "floods += regions_by_country.flood_poly.str.cat(sep=\", \")\n",
        "floods += ']}'"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "0yR7ouKBSr1e"
      },
      "outputs": [],
      "source": [
        "# @title Extract the features from road_geometry\n",
        "highways = '{\"type\": \"FeatureCollection\", \"features\":['\n",
        "highways += regions_by_country.road_geometry.str.cat(sep=\", \")\n",
        "highways += ']}'"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "9zQH6LjRV8Y2"
      },
      "outputs": [],
      "source": [
        "# @title Display Flooded Areas and Highways using geemap\n",
        "import geemap\n",
        "from ipyleaflet import GeoJSON\n",
        "import json\n",
        "\n",
        "styling = {\"color\": \"red\", \"fillcolor\": \"red\"}\n",
        "\n",
        "flooded_areas = GeoJSON(\n",
        "    data=json.loads(floods),\n",
        "    name='Flooded areas'\n",
        ")\n",
        "\n",
        "flooded_highways = GeoJSON(\n",
        "    data=json.loads(highways),\n",
        "    name='Flooded roads',\n",
        "    style=styling\n",
        ")\n",
        "Map=geemap.Map()\n",
        "Map.setOptions(mapTypeId = 'HYBRID', styles = {}, types = [])\n",
        "Map.centerObject(aoi, 12)\n",
        "Map.add_layer(flooded_areas)\n",
        "Map.add_layer(flooded_highways)\n",
        "Map"
      ]
    }
  ],
  "metadata": {
    "colab": {
      "provenance": []
    },
    "kernelspec": {
      "display_name": "Python 3",
      "name": "python3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}